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PREFACE 


The David W. Taylor Lectures were initiated as a living memorial to 
our founder, in recognition of his many contributions to the science of 
naval architecture and naval hydromechanics. His systematic investiga- 
tion of resistance of ship hulls is universally known and used, but of 
equal importance was his use of hydrodynamic theory to solve practical 
problems. Many of the experimental techniques which he pioneered are 
still in use today (for example, the use of a spherical pitot tube for 
exploring the structure of a wake field). The system of mathematical 
lines developed by Taylor was used to develop many designs for the Navy 
long before the computer was invented. And perhaps most important of 
all, he established a tradition of applied scientific research at the 
"Model Basin" which has been carefully nurtured through the decades, and 
which we treasure and protect today. 

These lectures were conceived to support and strengthen this 
tradition. We will invite eminent scientists in fields closely related 
to the Center's work to spend a few weeks with us, to consult with and 
advise our working staff, and to give lectures on subjects of current 
interest. 

It is most fitting that Professor Reinier Timman, mathematician and 
philosopher, initiate this series. He has long been a friend and on 
several occasions has used the Center for a retreat, to his benefit and 
ours. He has inspired and advised our staff and cooperated in our work. 
His students at Delft have made leading contributions to the development 
of modern naval hydrodynamics. Professor Timman's belief that mathe- 
matics can contribute powerfully to our technology is much in the David 
Taylor tradition. We are honored that he agreed to give the first in 


this David W. Taylor Lecture Series. 


W. E. CUMMINS 
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FOREWORD 


It is great honor to me to be invited to give the first 
in the series of David W. Taylor Lectures. My associations 
with the Model Basin date from a long time ago, and a visit 
to the United States is for me not a real visit unless I 
have the opportunity to taste once more the stimulating 
atmosphere which not only gives the Model Basin an out 
standing place in hydrodynamical research but also acts as 
a breeding ground where nearly all outstanding people in the 
field passed an essential period in their lives. So I am 
extremely grateful to have been given the opportunity once 
more to spend some time at this most interesting place and 
to participate in its work. I wish to express my gratitude 
to Justin McCarthy who originated the idea of the lectures 
and to all other friends who made this period a success. 

In particular, I am pleased that Dr. Langan, whom I used to 
know as a promising undergraduate student, did a fine job 


in editing the lectures. 


R. TIMMAN 


. 
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ABSTRACT 
The lectures present an introduction to modern control 
theory. Calculus of variations is used to study the problem 
of determining the optimal control for a deterministic sys- 
tem without constraints and for one with constraints. The 
method of dynamic programming is also used to solve the 
unconstrained control problem. Stochastic systems are intro- 


duced, and the Kalman-Bucy filter is derived. 


INTRODUCTION 

Optimal control theory is involved with the great human effort to 
control or influence processes of one type or another. The objectives 
and criteria for the performance of a physical system may be diffused or 
defy tractable analysis in many situations, but the basic concepts on 
which to proceed have been established in control theory. One first 
considers a system and a process through which the state of the system 
is changing in time; in other words, some action or motion of the system 
takes place in time. This behavior of the system is described by a set 
of time-dependent variables x (t) = (x) > a 3S x) which are called 
the state variables. In addition to the state of the system, one also 
considers controls by which the process in question can be influenced. 


These controls are represented by a set of variables u (t) = 


(uy Ge) eis un (t)) which are called the control variables. 
At a certain instant in time, say to» the state of the system is 
known to be x,. If an analysis of the system is to be performed, a sys- 


0 
tem of equations must be specified which predict the state for t > ty 


and for a given control function u. These equations are called the 
dynamic equations for the system; they may take the form of an ordinary 


differential equation 


or a difference equation 


x =f 

nt+1 re mn? a) 

They might even take the form of an integro-differential-difference 
equation or a time delay equation, but they cannot take on a form such 
that, the solution at some time ty is dependent on the solution in the 
future, t > t,. The dynamic equations must reflect this principle of 


al 
nonanticipation. One does not violate this principle by choosing a 
control in anticipation of the future and thus influencing the future 
state of the system based on estimated future information; in fact, the 
choice of such a control is actually based on the history of the state 
of the system available at the time of the choice. 

If no further specification of system performance is given, every 
control function which yielded a physical realizable state of the system 
for t > to would be a solution to the control problem. One can have a 
meaningful control problem only if there is a desired objective, a goal 
to be achieved by the process. Moreover, it is not sufficient merely to 
have a goal; there must be a control by which this goal can be achieved. 
This (controll Could be the case of no controll, £(t, x, 1) = (es =)5 
however, it must exist. Since it is not the purpose of these notes to 
delve into all the mathematical problems, it will be assumed that there 
exists at least one control by which the objective can be achieved. It 
will further be assumed that any control function used in the sequel 
yields a unique state function x (t) with x (to) = X35 the state func- 
tion is obtained by solving the dynamic equations. 

In general, there are a number of controls which could yield the 
desired system state. From among this set of possible controls, one 
would like to choose the "best" control with respect to some performance 
criterion. For example, one would like to choose the control so that 
the process is carried out with a minimum cost in fuel, or time, or 
money. In the sequel, it is assumed that the performance criterion can 
be expressed in terms of a cost function; furthermore, it is assumed 


that the cost function is additive with respect to the contribution from 


each time interval. An example of such a cost function is 


T 
G(x, T) + i F(o, x, u) do 


fo 


where Xp 


the system through the function G and on the intermediate states and the 


= x(T). This cost function is dependent on the final state of 


control function through the function F. The additive property of the 
control function with respect to the intermediate times is represented 
by the integral. By an optimal control is meant that control which 
minimizes the cost function; it is this function which is the desired 
result of optimal control theory. 

Any process that is being controlled is subject to unpredicted 
disturbances, and these can make a significant difference in the choice 
of a control function. Suppose the dynamic equations of a system is 


given by the differential equation 


where p(t) represents a disturbance. The behavior of the system in 
response to the two different controls (uy, = - x) and (uy =- ape) does 
not differ if there is no disturbance (p = 0); however, if a disturbance 


is present, the response is significantly different. If x, = 1, the 


0 
response to the first control is given by 


whereas the response to the second control is 


ie 
san= ang + | p(o) do 
« 0 


Such differences could conceivably result in a different choice for an 


optimal control. 


In analyzing systems and their control, one must find a way to 
represent the unpredictable disturbances. Such disturbances cannot be 
modeled by analytic functions since the value of an analytic function at 
any point is predictable from its value on an arbitrary short interval. 
One answer to modeling these disturbances is to describe them as stochastic 
sEoeaeses.” The theory of such processes was developed to model the 
fluctuation observed in physical systems. Wiener processes or the 
Brownian motion process are of particular interest to the stochastic 
control problem; many of the disturbances that affect a control system 
can be modeled by processes generated from Wiener processes. A Wiener 
process is a stochastic process in which the statistical properties over 
the interval (t, ttt) are the same as those over the interval (s, s+t); 
moreover, the behavior of the process is independent over time intervals 
which do not overlap, and there is no trend in the behavior. 

Once the stochastic disturbances have been introduced into the 
control theory, the problem is no longer deterministic. The state 
variables and control variables are no longer predictable but must be 
described by their statistical properties. Kalman and Bue provide a 
solution to the stochastic control problem for nonstationary linear 
systems. Their solution consists of using an optimal filter to estimate 
from the observed system performance the state of the system in terms of 
the conditional mean; the estimated state is fed back to the control 
signal through linear feedback. The linear feedback is determined by 
solving a deterministic control problem; the filter depends on the 
disturbances and on the system dynamics, but it is independent of the 
cost. Although the nonlinear stochastic control problem or its equiva- 
lent, the nonlinear filter problem, has not been solved, some headway 
has been made by Bucy and Jeaeping this lecture considers only the 


linear problem. 


Aslerem, K. J., "Introduction to Stochastic Control Theory," Academic 
Press, Inc., New York (1970). 


alle, R. E. and R. S. Bucy, "New Results in Linear Filtering and 
Prediction Theory," Journal of Basic Engineering Series D, American 
Society of Mechanical Engineers, Vol. 83, pp. 95--108 (1961). 


“ney. R. S. and P. D. Joseph, "Filtering for Stochastic Processes with 
Applications to Guidance," Interscience Publishers, Inc., New York (1968). 
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As an example of a control problem, consider a ship moving through 
a current of water; the ship is a system undergoing a change in state. 
In this example, the state is the position (x, y) of the ship. The 
parameters which control the motion of the ship are the power, which 
determines the velocity relative to the water, and the steerage angle, 
which controls the heading angle 8. In this simplification of the 


system, the dynamic equations are: 


be 
il] 


VN COS Gita Gxemmnyp) 


V sin 6 + v(x, y) 


S| 
Il 


where u and v are the velocity of the current in the x- and y-directions, 
respectively. The goal might be to go from point A to point B. If it 

is desired to reach B in the shortest possible time, the cost function 
would be the accumulated time; if it is desired to reach B with the 
minimum expenditure in fuel, the cost function would give the expended 
fuel in terms of x, y, V, and 8. A more complicated cost function would 
result if it is desired to reach B in the least time with a reasonable 
expenditure of fuel. Both the power and steerage angle could be subject 
to unpredictable perturbations; there could also be a stochastic pertur- 
bation of the current. 

This lecture on control theory first treats a deterministic optimal 
control problem with no constraints on the controls. It is first solved 
by transforming the problem into a boundary-value problem for an ordin- 
ary differential equation, the so called indirect approach; it is then 
solved by the direct method developed by Bellman, the method of dynamic 
a reaptemmeenee The big contribution of modern control theory to the de- 
terministic control problem has been the extensions to controls with 
constraints, and a discussion of constrained controls constitutes an- 
other major topic of the lecture. Still another important area is the 


development of the theory of stochastic processes necessary in the 


“Heilman. R. E. and S. E. Dreyfus, "Applied Dynamic Programming," 
Princeton University Press, Princeton, N.J. (1962). 


treatment of stochastic controls. Finally, the theory of Kalman-Bucy 
filters is given and their solution to the stochastic control problem is 


presented for linear systems. 


THE OPTIMAL CONTROL PROBLEM 
In these lectures the simplest optimal control problem considered 
is that of a state variable x(t) and a control variable u(t) defined on 
an interval O<t<T. The process being controlled is described by the 


dynamic equations 
BE) = (Ey Bo W) Gigi) 
with 


x(0) = x GR) 


0 
The vector f is twice continuously differentiable with respect to x and 
Lipschitz continuous with respect to u; this latter condition means 
simply that there is a constant L such that for every pair of control 


vectors u and v 
IEG, & w) = 2, => wD) <aloe x (1.3) 


For each control vector u, these conditions imply that the state 
vector x, which is obtained from solving (1.1) and which also satisfies 
the initial condition (1.2), exists and is unique. Moreover, from among 
the set of control vectors, it is assumed that there is a unique control 


u which minimizes the cost function C The cost function is defined by 


T 
the following: 


T 
CS lw) =sCl Ga) ao | F(o, x, u) do (1.4) 
qt =r) A xu 


The functions F is twice continuously differentiable with respect to x 
and Lipschitz continuous with respect to u; G represents the cost at the 


terminal point x(T) = x it is twice continuously differentiable with 


Tt? 
respect to Xn 
Suppose that v is an optimal control vector, and consider a slight 


deviation Ou of this control vector. If 
u(t) =v + du 
u(t) is also a control vector, as can be seen from an application of the 


theory of ordinary differential equations. If z is the state vector 


associated with the control v, the new control u yields a new state 


vector x given by 
x(t) = 2 + 6x 


where 6x is an unknown. Moreover, since v minimizes the cost function, 


the new cost function is greater; 


al 
1 (Of: Cs (OM a Cx @ cer ID) 
B =u 7 


IE 
> { D@, 2 w) d@ + G(Z pT) (@>)) 
0 


Since the old state vector satisfies 


Z 


tl 
Fh 
= 
in 
. 
N 
. 
<q 
Vv 


and the new one satisfies 


| 
iT] 
Fh 
oo 
ct 
x 
(ory 
4 


then 


f(t, z + dx, v + du) 


Now by assumption f is twice continuously differentiable with respect to 


x; hence 
ox = f(t, z+ ox, v + du) - E(t, z, v) 
= £, dx + £(t, 2, v + du) - £(t, z, v) + 0(|dx|*) (2.6) 

It is not necessary that du be uniformly small; indeed, in problems 
involving bang-bang controls, this is not at all true. However, there 
can be deviations du of order one only if their duration is short. It 
can be proved that if du satisfies the condition 

Ge) 


T 
{ |Su(o)| do <e 
0 


Since by assumption, f 


then the deviations 6x(t) are also of order €. 


is Lipschitz continuous with respect to u, 
= 0(du) 


l(E5 25 B) > Ges Bs w)) <ulo = a 


Moreover, it follows from Equation (1.6) that to the same order of 


approximation 
(1.8) 


Ses se Os ae (ES 5 Wy) = (65 5 WY) 


or in abbreviated form 
(1.9) 


Ox = £ Oxi E(WW)NS FW) 


This equation is a linear differential equation for 6x, and there 
é ; 7 ‘ F 5 ; 
are standard ways for solving linear differential equations. One first 


considers the linear homogeneous equation 
yay (1.10) 
where in our case y represents the vector 6x and A the matrix Se Let 
yi(t) = Cis Der eae an a (ene 
be the solution of Equation (1.10) with Bee (Ee T) = OEP the Kronecker 


delta; moreover, let @(t, tT) be the matrix whose column vectors are the 


vectors ee O(t, T) = Bes T). The matrix ®(t, T) is called the 


transport matrix or fundamental matrix for the differential Equation 


(1.10). From (1.10) it follows that as a function of t 


= Ge. ©) Shite, © (1.11) 


and by its definition 
OE, t) = I Gea) 


where I is the unit matrix. The solution y(t) is given in terms of its 


value at t = T by 


SAGE) = OCES a) evar) Giesiss)) 


oeeddiae tod E. A. and N. Levinson, "Theory of Ordinary Differential 
Equations," McGraw-Hill Book Company, Inc., New York (1955). 


Hence 


USGS) S08) = CE a) a7 (Ce) 
= Ol GE) se) Crement=)) mya Ct) 
or if y # 0, 
TL =O (5p) eOlGee se) (1.14) 


Differentiating with respect to t yields 


o® 


Cie es 


(t, 7) @(t, t) + O(t, 1) Z O(t, t) 


A(t, 1) O(t, t) + O(t, 1) $ act, ¢) 


A(t) + (t, 1) £ (2, ¢) 
It can be shown that 9(t, Tt) has an inverse and that this inverse is 


O(t, t); consequently 


4 OG, Boo OE wD Ate) = = OG, &) AG) 


that is, ®(t, T) as a function of T satisfies 


oC, T) = - O(t, T) A(t) (1.15) 


Although (1.15) will be used subsequently, of immediate interest is the 


solution to the inhomogeneous linear equation 


y = Ay te (ee) (1.16) 


10 


with y(t) = 0; the solution is given by 


t 
y(t) = { O(t, 0) g(o) do (il, 1177) 
1E 
which can be verified by substitution into (1.16). For the control 
problem, (1.17) has two consequences: it can be used in conjunction 
with (1.6) to obtain an estimate for the order of magnitude of 6x and 
it can be used to solve (1.9). 


In the first case, 


e t 


| dx | i [o(t, o)|| £(@) - £~)| do + { 0(6x7) do 
0 0 


[A 


ie t 


M { | £(u) - £()|do + \ 0(8x7) do 
0 0 


[A 


where M is a bound for ®. From (1.3) 


t 


t 
| Sx| LM i | Su| do + { 0 (6x2) do 
0 0 


JA 


t 
IMe + | 0(éx") do 
0 


[A 


By iteration 
t 
|ox| < LMe + { 0(e7) do = O(c) 
0 


The second case is of more interest, of course, for it gives an 


approximation of 6x good to the second order in €, namely, 
t 


6x = il (Es (6) ea) FG) lindo (1.18) 
0 


11 


where ® is defined by 


SP (t,t) = £ (te) O(e, 1) (1.19) 


Now consider the difference in the values of the cost function; by 


GE) 
(E 


il D@s, > W)) = BW], 45 W) do + C(xns T) - C(Za> 2) = © 
0 


Hence, from the assumptions on F and G, 
TE 


{ [P_ OS se EG) = B@)]) cde + E_Ox(t) > © 
J oo ry re a 


By (1.18), 
T 


iD 
{ ERA) i OG, CS) Gao) > £@©) pido Cae) riGiGe)) ieee 
0 0 


aE 
+6 { 0(T, o) (E(u(o)) - £(v(c))) do > 0 (1.20) 
0 


If the order of integration in the double integral is changed, 


it T 
{ Ba Gp) i OG nO) LEG @) eet Gio) adodyn 
Oi 0 
T T 
= i i F(t) O(t, OC) di Guo) eae Gn) sido Glo ZL) 
0 o 


The vector function p. is defined by 
a 


pute) & = | a) OG, =) tr = CG) O@, &) G22) 
t= 


Recall that one of the properties of ® was (1.15) 


oe (t, t) = - O(n, t) £ (ct) 


12 


Then 


Te 
ie = F(t) O(t, t) - i F(t) “ O€e, ie) ae = G(T) a Ol Ghee) 
t 
TE 
= F(t) + J F(t) O(T, t) f(t) dt + G(T) o(T, t) f(t) 
ap 
= F(t) - J F(t) OG) dit — G(T) O(T, 0 f(t) 
B = sj £ G23) 


with p(t) = G(T). In terms of a. (1.20) becomes 


ili 


| [- p'(c) (E(u) - £(v)) + Fu) - F(w)] do > 0 
0 


or 
il T 
| l-2e@ rp EMI ole ?@ rp E@ de2O0 Gers) 
0 
Since du is an arbitrary deviation satisfying only (1.7), it can be 


chosen such that u = v everywhere except on some arbitrary interval; as 


a consequence, the inequality in (1.24) must hold for the integrand: 


- Fv) +p) £@) > - F@) tp EW) 
Define 
H(t, u) = - F(u) +p’ £(u) (1.25) 
Then H satisfies 
H(t, v) > H(t, u) (1.26) 


13) 


for v, an optimal control. This is the Pontryagin maximal principle 

3 ; T ; : 
which states that for given values of p and x at time t, the optimal 
control v(t) is the control function for which the Hamiltonian H(t, u) 


is a maximum. 


If the control functions are sufficiently smooth, the optimal 


control is that control for which 


je hs Sp eps (1.27) 
u u 


It is assumed that f is differentiable with respect to u; prior to this 
equation, f need only be Lipschitz continuous with respect to u. This 
equation is a system of m equations which could be solved for the 

m control functions (uj5+-+-,u_) in terms of the state variables 
(X)5-+-+-5%)) and the new variables (@iocoesd yc Consequently, the 
optimal control problem has been reduced to a two-point, boundary-value 


problem for an ordinary differential equation: 


oe S CE, Xs u) 


Die hreDe oe 
T 
=- + 
0 i Dp i 
or 
5 oH 
x = = 
op- 
oH 
.ea56 e 16 2S 
p = ( ) 
_ OH 
ons du 


14 


where 


x(0) = Xo 


pe) = G(T) (1.29) 


There are just enough conditions to determine x, pi, and u. 

The function H contains the variables x, rae u, and of course t. 
Using (1.26) to eliminate u, (1.28) can be expressed in terms of the set 
of dual variables x and p' = Ds where the prime denotes transpose of 
the vector; the resulting system is the familiar canonical form of 


classical mechanics. 


|e 
| 
| 


p=-— isa) 


The boundary conditions are stated in terms of Xo» Xp» and T; for 
instance, both Xp and T might be fixed, or either one might vary while 
the other is fixed. No boundary conditions are specified directly in 
terms of p; the boundary conditions on p are obtained indirectly by 
substitution into (1.29). Equation (1.29) does, however, contain a 
sufficient set of conditions to pose a two-point, boundary-value problem 
for (1.30). 

Another form that the boundary condition at t = T might assume is 


for x, 


T and T to satisfy an end condition of the form 


M(x, T) = 0 (1.31) 
where M is a twice continuously differentiable vector function of both 


its arguments. In this case the method of Lagrange multipliers will be 


used to transform the optimal control problem into a corresponding 


WS) 


two-point, boundary-value problem. The vector q is introduced here as a 
Lagrange variable. Now the problem of minimizing the cost function (14) 


is replaced by the problem of finding the unconstrained minimum of 


iP 
Ce @ = { F(o, x, u) + q'(x - £(0, x, u) do 
0 


+ IMGs T) + C(x LD) AGES 29) 
The boundary condition (1.31) has been inserted into the cost function 
by means of the Lagrange multiplier uw. Suppose v(t) is the control 
which minimizes Gao For a variation du to the control v, let x(t) 
denote the new state variable, and let t = T+AT be the time at the new 
terminal point. The main difference from the previous argument in this 
section is that the terminal time is T+AT rather than AT. The new cost 
is given by 

T+AT 
Cw = N@_ Bsr Os, War Ow) GY]! sw OR = (> Bw Ok Wo Su) do 


+ G(x(T + AT), T+ AT) + u'M(x(T + AT), T + AT) 


Hence the increase in cost C2 - {2 is given as: 


ue 


DOE Bs a eel 
J [F dx + Fou + q'6x q'f ox q £ oul do 


ca = OD) 


37 y 
+ G Ax, + GAT + u"(M Ax, + MAT) 


cee 


T+AT 
| (On xe) tach = igi (O75) x.) tu) ido 
2 


where x = z + 6x 


EN ree T 
q' Ace dt = q'(T) 6x(T) - q' 6x do 
Once a iF aaa a 
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Hence 


lee 


QO < 4 WES = q' - g'f] Oxait [F - gine du} do 


+ ' f x i] 
q' (1) 6x(T) + G Ax, + GAT + u"(M Ax, + MpAT) 
T+AT . 
fs | F(o, x, u) + q'x - g'f(o, x, u) - F(T, & uD) 
a 


= GsK(GE) se GTEC, xx 


Xp» up) do 


ar (USE. Xp uy) + q'Zn aF @iYOxe(Gb) = @ VEC Xr up) | AT (1.33) 
The integral from T to T+AT is a second order contribution which goes to 
zero faster than the other terms as AT > 0. 
In order to determine AX» consider the solutions of the differ- 
ential Equation (1.1), which have the initial value x,. These solutions 


=0 
satisfy the integral equation 


t 
x(t) = X> + | £(0, x, u) do 
0 
Then he 
Noon (Ga) zi) | E@; = w) do 
ae bne a) T 


= 6x(T) + £ AT + 0(e*) 


where the geometry of the proof is illustrated in Figure l. 


A T+AT 
Figure 1 -- Geometry of the Proof 
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To within second order 


6x(T) = Ax, - £ AT (1.34) 


Within this order of approximation, (1.33) reduces to the following: 


ae [ Sian! q coll Ope AP bs qf du do 
+ [q' ' 
fey CE) sh Cae WE Aas, 
AES) ' 
qe [l= @" Cit) ze sp CG. ar iu My + F(T, Xp uy) | AT 
If gq is now determined so that the coefficient 


of 6x vanishes, 


This is the same differential Equation (1.23) that p satisfied; our 


Lagrange multiplier can then be identified with p 


a> (1.35) 


Moreover, since the relationship must hold independent of du, 


= ' oo 
FY pif, 0 (1. 36) 


Since there are no longer restrictions on Ax and AT, 


(is Ss) 


Introducing the Hamiltonian (1.25) yields 
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a= 3p 

5 oH 

‘pe as Dee (1.28) 
_ 3H 

Oe du 


The initial condition x(0) = x, together with the terminal conditions 


—0 
M(x, T) = 0 
DCE) ES aaa (uy Mit Ce) Giese) 
Gs mC) = G,  woM. 


provides a sufficient number of conditions to determine x, p, u, and T. 
The last two equations in the system (1.38) are obtained from (1.37). 

The problems of optimal control theory generally reduce to a two- 
point, boundary-value problem for the system of ordinary differential 
equations (1.30). Bailey, Shampine, and Walemanc discuss methods for 
solving such two-point, boundary-value probiems. These problems are 
presently solved either by the shooting method or by solving a sequence 
of simpler boundary value problems whose solutions converge to a solu- 
tion of the given problem. In any case, very few of these problems can 
be solved without the use of electronic computers either digital or 
hybrid. 

The shooting method is the easier, when it works. It consists of 
supplementing the conditions at one end with a sufficient number of 


assumed conditions to yield an initial value problem. The initial value 


yale, P. B., et al., "Nonlinear Two-Point, Boundary-Value Problems," 
Academic Press, Inc., New York (1968). 
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problem is solved; the solution is substituted into the boundary con- 
ditions at the other end. If these conditions are satisfied, the solu- 
tion to the initial value problem is the desired solution to the two- 
point, boundary-value problem; otherwise, a new set of assumptions is 
made based on the discrepancy between the actual boundary values and the 
calculated values. Hopefully, as one continues this iteration process, 
the solutions to the initial value problem converge to a solution of the 
two-point, boundary-value problem. The shooting method may not con- 
verge, or it can be unstable, that is, a small variation in the initial 
conditions results in a large variation in the solution. If the initial 
problem is unstable, a small error, such as roundoff on a computer, 
could cause subsequently computed values at another point to be meaningless. 
Before proceeding to the direct method for solving the optimal 
control problem, take a second look at the Hamiltonian H and the func- 
tions nice Suppose that the terminal cost G is identically zero; the 


cost function is then 


T 
C(u) = | F(a, x, u) do 
0 


Further, assume that every point in an open neighborhood N of an optimal 
trajectory z(t) can be joined to the initial point (0, Xo) by a trajec-— 
tory x(t) resulting from an optimal control. This assumption makes the 


minimal cost J a function of the terminal point (T, X) in N. 
a 
J (p> T) = Min i (Gs x5 Ww) ao (1.39) 
0 


It is assumed that J is twice continuously differentiable. Then, 


J (Sp + AX T+AT) = J (&n> T) + J Ax + JpAT (1.40) 
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By the definition of J, there is a control u + du together with a 


trajectory x + 6x such that 


als 
J (Sp + AXn> T+AT) = { R(@5 soar Oey Woe Su) do 
(1.41) 
where u is the control such that 
T 
Js T) = | WW, 35 wy) do (1.42) 


0 


From (1.41) and (1.42) 
J (Xp + AX» T+AT) - J (ips T) 


T 


= { I(@, a OX> War Ow) = F(, x, wu) do 
0 
T+AT 
+ N@, + OX, Ww OW) do 
T 
Ie 


{ (F 6x + F du) dt + FAT 
0 x — nm — 


Now from (1.23), the equation for p is 


Hence 
J (Kp + Axon» TFL) = J (Xs Ty) 
aE 
fa [tq + pte) ox + F bul do + FAT 
0 SK Ww = 


T 
peGL) soscGh) ace J p'(£, 6x - OkS)ch F ou do + FAT 


ZAM 


From (1.9) 


Hence 


_— = ' — 
J Sp + AX T+AT) J Gn> sl) p'(T) (Ax, f£AT) 


T 
lane 
; @, = pie) Sa ce <p mn 


where use has been made of (1.34). But by (1.27), FY pot = 0; so, 


JX + Axas T+AT) - J(x T) = p' Ax, ae (OB = jy (ie)ie) Aa 


T°? 


p'(T) AXn - HAT 


ad ite, b gin 
where the last equality results from (1.40). This gives 


J =p (1.43) 


and 


Jp 2-H (1.44) 


In the space of variables (Xp T), the vector p' is the gradient of 
the function J; it is normal to the surfaces of constant J; H is the 
Hamiltonian of the function J. This sheds new light on the maximal 
principle. Along an optimal trajectory, the change in cost J over a 


given time step AT is a minimum, that is, H is a maximum. 
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OPTIMAL 
TRAJECTORY 
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| 
|<— AT—>| 


| 
| 
| 
| 


Figure 2 -- Constant Cost Fronts 
These arguments hold only if the terminal cost is zero; G = 0. 


RELATION TO DYNAMIC PROGRAMMING 

The partial differential equations (1.43) and (1.44) can be ob- 
tained by the method of dynamic programming. This method is based on 
the Bellman principle of optimalicy.! According to the Bellman prin- 
ciple, an optimal control policy has the property that, regardless of 
the initial state or initial decision, the remaining decisions must 
constitute an optimal control policy with regard to the state which 
results from the first decision. 


In terms of the cost function 


A 


C(u) = i F(G, x, u) do 
0) 


the Bellman principle takes the form. 


The cost C(u) is a minimum along a curve x defined on [0, TE] suie abie 


is a minimum along each later part of the curve, that is, if 


Varese. S. E., "Dynamic Programming and the Calculus of Variation," 
Academic Press, Inc., New York (1965). 
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oT 


il F(o, x, u) do 
t 


is a minimum along the curve x on the interval [t, T] for all te[0, T]. 
The integral is dependent on the end point (t, x(t)). If one 


defines 
ab 


J(x, t) = min | WG, 385 wy) do @m) 
u t 


for all admisible controls u, then 


T 
Je, te) S malin Wee, se u) 6t} + min { F do 
u UW EOE 


or 


JiGs) tt) = min THe, x. uw) Ot + JiGg + Ox, tt Ob TCR 
u 


This equation forms the basis of the direct methods for solving control 
problems, described by Teerines! Lamson extended the direct methods to 
constrained problems. 

If it is assumed that J has partial derivatives, the differential 
equations (1.43) and (1.44) can be obtained from (2.2). Hence, the 
boundary value problem for the optimal control is obtained. If the 
partial derivatives of J exist, the right-hand side of (2.2) can be 
expanded in a Taylor series: 


J(x, t) = min {Fot + J(x, t) + J Gx, i) Ox ou J x, ©) Orel 
u 


(2.3) 


Sparse, R. E., "State Increment Dynamic Programming," American 
Elsevier Publishing Company, Inc., New York (1968). 
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From the differential equation (1.1) 


Hence 


0=min {F+J£+ I} &t 
t iS ic 


Since dt > 0, 


On= mae {F + Jf + Jt CES) 


In order to find the minimum of the term in brackets, it is differ- 
entiated with respect to u and the result is set equal to zero. This is 
a necessary, but not sufficient condition; however, if one assumes a 


minimum, it serves the purpose. 


oe alamo) 0) (2.5) 
u x— 
Biya (Zi) 
we wy iad, = © @G)) 
x t 
Hence 
H =0 
u 
(1.43) 
JF ao BS = (i a Jf) 
From: (1.25) 
AL 
J, =p (1.42) 
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There is a difference between the definition of H here and its 
definition in the previous section. This is only an apparent difference 
in the sign of F, which occurs because the lower limit of the integral 
is used in the definition of J here rather than the upper limit as used 
earlier. Otherwise there is complete agreement with the results of the 


indirect method. 


CONSTRAINTS ON THE CONTROL AND STATE VARIABLES 
In most applications, the control or the state variables cannot be 
chosen arbitrarily but are subject to constraints. In the problem of a 
ship moving in a current, ship speed is limited by the maximum power 
available. The constraints can generally be expressed in terms of 


inequalities of the form 
g(x, u) < 0 Ca» 


where the vector inequality simply means that the components satisfy the 
inequality. The number of components in the vector ¢ is the number of 
constraints on the system. The analysis does not depend on whether both 
x and u occur implicitly in the inequality; one can have constraints on 
the controls and not on the state of the system or vice versa without 
affecting the analysis. 

In this presentation, the variables in the optimal control problem 
with constraints are the state variable x(t) and the control variable 
u(t) defined on an interval 0 < t < T. The process being controlled is 


described by the dynamic equation (1.1): 
a(S ECE, ze WW) 


with initial condition x = x5 the state and control variables are 
constrained by the inequality (3.1). For simplicity, the terminal cost 
is taken as zero, G = 0, and the cost function is given by the equation: 
T 
C,(@) = J (5 285 wy) Glo Ge) 
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The vector f and the cost function F are twice continuously differentiable 
with respect to x and continuously differentiable with respect to u. 
The Lagrange multipliers will be used here to reduce this problem 
to a two-point, boundary-value problem. As in (1.32), the differ- 
ential equation is introduced into the cost function by means of a 


Lagrange multiplier p. 


T 
I@, 35 BW) py Ce = 2@5 &5 W)) ae 
0 


os (u) 


which yields the variational equation 


ap 
i 2 te 1! Ow ar)? Oe = jo)! 
0 gia tae aa ay 


Pe) 
So. £ ox Die Su] do 
T 
= ' =a iT 
= p'(T)éx(T) + J [@, - B' - p'£,) ox 
as ' 
+ = pie) & do 2 O (3.3) 


The differential in the cost is greater than or equal to zero since 
it is assumed that the variation du is around an optimal control, a 
control which minimizes the cost. 

Because of the constraint (3.1), the vector du is not free. 
For instance, suppose that for t between ty and tos the trajectory z(t) 
due to the optimal control v(t) is along the boundary of the allow- 


able region; see Figure 3. One cannot freely choose the variation Su 


in the control vector fort, <t<t 


1 and still expect to remain in the 


2 
allowable region R. 
For the optimal trajectory z and control v, there are at most a 


finite number of intervals t, < t < t, + 1 such that equality holds for 


k k 
any of the equations in (3.1).* On such an interval, the conditions 


(3.1) can be split into two sets 


*The proof of the statement is topological and beyond the scope of 
these notes. 
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—_ z(t) 


Seocooe x(t) =2zat 6x 
Figure 3 -- Constrained Variables 
(Wp XQ) = Ww 
and 
o, (2, w) < O (3.4) 


where > = (b> d5)- 


Consider a new vector defined by 
o(x, u) + px, u) = 0 


The vector ) is called a defect vector. Along the optimal trajectory, 
the vector can also be split into two component vectors, vy and Wo> 
which correspond to the component vectors of ). The component vectors 


of p~ also change from interval to interval. Along a given interval 


Les Sa 
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(3.5) 


Since , v) is zero on this interval, either z, v, or both are 
on the boundary of their allowable range. From previous arguments, it 
is known that one cannot freely choose du. Only those values of du are 


allowed which satisfy 
o, (2 + Sx, v + du) < 0 
or by (3.4) 


go, (2 + Ox, v + du) - $ (2, vy) < 0 
On the other hand, for a neighboring trajectory to z 
ott OS War OH a Ws Ob = O 


on [t,> thayl- Simca.) 


WE a OX Wa Ot) = Os, W) = Ol = © (3.6) 


In order that the above inequality and (3.6) hold, 
sy, 2 0 (3.7) 
Moreover, provided the variations are sufficiently small, I) is free. 


If » is twice continuously differentiable, then it follows from 


(3.6) that 


>, Sx +o, Sut Sp =O (3.8) 
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Set du = (Su, > du.) and consider the first N, equations in (3.8), 


OF) 
Xo, = dim (o,)- 


Sy, OF Su iu, Sup? SU ae 


If the square matrix is not singular, its inverse exists, and 
ile 8 ’ Vf 
1 


Co a iu, OB TY ise CS 


The vectors Su, and 6x are free; the vector ob satisfies (3.7). If 


the matrix Oa is singular, the first N constraints were de- 


>) 
pendent; eliminate the dependent constraints and start again. 
The contribution to the cost differential (3.3) from the 


IMBeIVEL (ES ie Se is the following integral: 


k k+1 
(ee 
= oe (ae ie 1 at = 1 
te z {LF P P ee a Y ee DE, Y Wis 6x 
t aL 
k 
=a es = wy! 6 
ee Pa, "a pif) ¥ 9, 1 %% 
a5 ' 
=. p'é ) ¥ SH) do 
iL alk 
Define the vector Ay by 
Mes SOR Sie) Ay (3.10) 
% AT Cane 
Then 
oc 
we = Papen ' 5 
De = {LEY Pp pase POs call O. 
Te 
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The vector p can be determined so that the coefficient of 6x vanishes: 


sui + al U pple 
= Be Be. © Bee oo 
Since du, is free, the usual argument that Su, is zero everywhere except 


on a small interval yields 


Po op stl oe = © (3.12) 
2 2 


Now du can be chosen so that du = O for t < t, and for t,,, <t. 


In this case, the only contribution to the cost difference (3.3) is that 


due to I hence 


(ee 

= 1 

ORS AG oy do 
Fig 


By (Boo Op, > 0; so 


> 0 (3.13) 


Let the Hamiltonian be defined by 


== 1 ap jp! 2 > AY © (3.14) 
where X is defined by 
KN, © © she @, = © 
4 
Mp, =O WE w, < YW 
J °5 
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The differential system (1.28) also holds for this H, that is, 


=, Gel 
x dp 
__ i 
DR Sie tae. (1.28) 
oH 
OS 


One example of a constrained control problem is that of a forced 
harmonic oscillator in which the magnitude of the force is limited. In 
this problem, the force is the control and the process is one of chang- 
ing the velocity and displacement of the harmonic oscillator. It be- 
comes an optimal control problem if one is interested in finding the 
force or control which reduces the oscillator from a given velocity and 
displacement to zero velocity and displacement in minimum time. 

The equation of motion for the forced harmonic oscillator with a 


limited force is simply 


where |F| <M, a given constant. Set x = cz/M, T = wt, and u = F/M 
where w = VYc/M. In terms of these nondimensional variables, the non- 
dimensional form of the equation of motion is 


x Pe Sw (3.15) 


where the control function satisfies the inequality |u| < 1. This 


constraint can also be written in the form 


32 


GpiGn) SeaiCaye D)esy0 
(Eee) 
On) == Cur ys © 


The optimal control problem can be formulated in the phase plane. 
If (x, y) are the phase plane coordinates, the equation of motion (3.15) 


takes the form 
(Gree) 


Starting the oscillator at a given displacement with a given velocity is 
equivalent to assigning a given point (x, y) = (a, b) in the phase plane 
as an initial condition for (3.17). The rest state of the oscillator is 
represented in the phase plane by the point (0, 0), the point of zero 
displacement and velocity. Hence, the optimal time control problem is 
one of finding a control u which minimizes the time between states 


(a, b) and (0, 0). In this problem, the cost is given by 


T 
C@) +P | dt (3.18) 


The cost function F(T, x, u) = l. 


Set p = (p, q). Then the Hamiltonian defined bye @rel)) ees 
H =- 1+ py + q(u - x) - A@u- 1) Gu t+ 1) (3.19) 


and, moreover, (1.28) takes the form 
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in ea ay Tey 
oH _ 

erat x ap {Ul 

=. oo 
6 5 = « 

eu ey SEE 
4 age Tae 
0 = Beg - A(u- 1) - A(u +1) (3.20) 


Suppose u is an optimal control which reduces the oscillator from 
the state (a, b) to the state (0, 0) in the minimal time T, and suppose 


|u| < 1 for the interval To St Siti Suppose q # O on T) Kt < T,. 


By (3.20), q - 2Au = 0; hence, X # O on (T,> T,)- A consequence of 


dX #0 is that » = 0; hence, if q # 0, it follows that Ju(t) | = lon (T,> 


In other words, one needs to look only for the optimal control among 


those controls for which Ju(t) | = 1. 
Now |u| Seeleetimplete satis —sehaels gunenceP mathe ms out onw ores (3-720) mers 
given as: 
x Ls A sim @ 2 a) 


A cos (t + a) 


< 
il 


> = 3) Sim Ce a Oy) 


) 


Bcos (tT+a 


QQ 
tl 


0 


eas (3.21) 


fe) 
i} 
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2) 


Since X > 0, it follows from the last of these equations that the sign 
of q is the same as the sign of u. Hence, if q changes from positive to 
negative, the optimal control must switch from +1 to -l. It switches 
from -l1 to +l if q changes from negative to positive. 

In a neighborhood of the origin, the optimal trajectory satisfies 


@ ED Sy =1 


Hence, its final segment is either on the circle of radius 1 about 
(-1, 0), or it is on the circle of radius 1 about (1, 0); see Figure 4. 
Suppose for the sake of argument that there is an € > O such that 
HG) se Foruh ie) San er themlasitysesment of ‘the joptimalsstragector 
is on the semicircle {(x + i)? + a = ils © < Wes 

Between (0, 0) and (-2, 0), the parameter T would change along this 
semicircle by the amount 1; hence, the sign of q must change somewhere 


on this semicircle. At the point S, where q changes sign, the sign of u 


must also change, and u switches cee, -l to 1. The optimal path continues 
backward on the circle of radius ry around (1, 0) until either (a, b) is 
reached or q changes sign. But q does not change sign until the point 

S5 is reached since the time between Sy and S5 I Wo | ANE So» the control 
would switch to -1l and the optimal trajectory would continue back on the 
circle of radius 5 around (-1, 0). This process is continued until the 
point (a, b) is reached. In the process, one switches control each time 


one of the following semicircles is intercepted: 
2 2 
be = Gme I) say = lo VS 05 MED SISZ5 e060 (@e22)) 


or 


ik @n oie Sr Sie SO BeO,,.4. GL25) 


The curve formed by these semicircles is called the switching curve; see 


Figure 5. 
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Figure 4 -- Optimal Trajectory 


Figure 5 -- Switching Curve 
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The optimal control and the resulting trajectory in the phase plane 
can now be obtained by reversing the above procedure. If (a, b) is 
above the switching curve, proceed with the control u=- 1. The 


optimal trajectory will be along the circle 


(x + 1)? + v = (at 1)? + be 


in the direction of that part of the switching curve which lies to the 
right of x = 0. For (a, b) on the switching curve, use u = - 1 if 
x < Q0oru=1if x > 0. If (a, b) lies below the switching curve, 
start with u = 1 and change to u = - 1 at the switching curve. 
Change the sign of u at each intersection with the switching curve. 
When u = 1, the optimal trajectory lies on a circle with center at 
(1, 0); when u=- 1, it is on a circle around (-l1, 0). 

Suppose only one switch in u is needed to reach the origin from 
(a, b). Because of the symmetry of the problem geometry in the phase 
plane, it is necessary to consider only those cases for which a = 1 


after the switch. The origin is then approached along the trajectory 


1 - cos (T - T) 


* 
Ul 


y = = sim (= Tt) (3.24) 


2 
which is on the semicircle {(x, y)|(x - Ly? yn —tollva Ole Met 
Ge be the time at which the switch occurs. The optimal trajectory 
for T Sous is given by 


== 1 Asim (G0) 


A cos (tT + a) (@EZ5)) 


< 
i] 
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where A and @ are constants defined by 


Asinag=adz+l1 


A cos @ = b 


By (3.24) and (3.25), the switching time must satisfy 


1 - cos (T - ue) 1 +A sin (T. + a) 


- sin (T - i) A cos (T. + a) 
Elimination of ue from these equations yields a relationship between the 


terminal time T and the initial point (a, b), namely, 
2 F 2 
(A cP iL sp Cos Tt)” a (> + otm WwW) = 4 (3.26) 


By definition, time fronts are the curves which connect initial 
points having the same terminal time T. Equation (3.26) can be used to 
determine the time fronts for T Soin LES Ole the eimesrrontetsms imple, 
the origin; if there are no switches in the control, the initial point 
is an endpoint of the curve connecting all initial points from which the 
origin is reached with one switch in time T. More than one switch would 
require T > 7. From (3.26), the time fronts for 0 < T < 7 are segments 
of the circle of radius 2 around the point (-l - cos T, - sin T); see 
Figure 6. It is the segment of the circle which lies above the switch- 
ing path. At the switching path, the time front is tangent to the 
vertical line x = constant for x > 0; at the opposite end, it is tangent 
to the switching curve. For T = 7, the time front is a circle of radius 


2 around the origin. 
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STOCHASTIC SYSTEMS 

Stochastic control theory was first applied in this country at the 
Massachusetts Institute of Technology during World War II to synthesize 
fire control systems. In the 1960's it was applied to space navigation, 
guidance, and orbit determination in such well-known missions as Ranger, 
Mariner, and Apollo. Applications of the filtering theory, aspects of 
control theory include submarine navigation, fire control, aircraft 
navigation, practical schemes for detection theory, and numerical in- 
tegration. There have also been industrial applications; one example 
involved the problem of basic weight control in the manufacture of 
paper. 

The filtering and prediction theory developed by Wiener and Kolmogorov 
forms the cornerstone of stochastic control theory. It provides an 
estimate of the signal or the state of a process on the basis of observa- 
tion of the signal additively corrupted by noise. Unfortunately, the 
Wiener-Kolmogorov theory cannot be applied extensively because it requires 
the solution of the Wiener-Hopf integral equation. It is difficult to 
obtain closed form solutions to this equation, and it is not an easy 
equation to solve numerically. 

Kalman and Buey- give a solution to the filtering problem under 
weaker assumptions than those of the original Wiener problem. Their 
solution makes it possible to solve prediction and filtering problems 
recursively and is ideally suited for digital computers. Basically, it 
can be viewed as an algorithm which, given the observation process, 
sequentially computes in real time the conditional distribution of the 
signal process. The estimated state of the process is given as the 
output of a linear dynamical system driven by the observations. One 
determines the coefficients for the dynamical system by solving an 
initial value problem for a differential equation. This differential 
equation is easier to solve than the Wiener-Hopf equation. 

Our attention here will be limited to linear systems with quadratic 


cost functions. In this case the solution of the optimal control 
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problem is given by the separation ihaoreemn = The solution consists of 
an optimal filter for estimating the state of the system from the ob- 
served data and a linear feedback of the estimated state of the system; 


see Figure 7. 
OBSERVED 


DATA 


PROCESS 


OPTIMAL 
FILTER 


CONTROL SIGNAL 


LINEAR 
FEEDBACK 


ESTIMATED STATE 


Figure 7 -- Stochastic Control System 


The optimal filter is the Kalman-Bucy filter, which will be dis- 
cussed in detail in the next section; the linear feedback is the same as 
would be obtained if the state of the system could be measured exactly 
and if there were no randum disturbances in the system. Thus, the 
linear feedback can be determined by solving a deterministic problem. 
Because of time limitations, we will not prove but merely accept the 
separation theorem. 

One objection to the use of stochastic control theory is that the 
process to which the theory is applied may not be random but merely 
irregular. For instance, the traffic flow on the Washington Beltway may 
not be truely random but it is certainly highly irregular. If I need to 
reach Dullis Airport from DINSRDC by 1 pm, it might take me 45 to 50 
minutes; but to reach the airport at 6 pm, I would have to allow 2 
hours. The reason for this variation in lead time is that there will be 
bumper-—to-bumper traffic on the Beltway during the rush hour and any 
accident brings this traffic to a halt. It is not the microscopic but 
the macroscopic properties of the traffic flow that govern our lead time 
estimate. The traffic flow could be analyzed as a stochastic process; 


such a model would be acceptable provided it predicted the macroscopic 
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properties of the flow. This is analogus to using linear models in the 
deterministic case. If the predictions agree with the experimental 
results, the linear theory is said to be good; if they do not, then the 
process is said to be nonlinear. In using a statistical model, one 
should recognize that it is only a model and not the actual process, and 
one should continually strive to determine the accuracy of his models. 

There are many reasons in favor of applying stochastic theory. The 
solution of the stochastic problem may be possible whereas the determin- 
istic theory may be hopelessly impossible. In many problems such as 
that of traffic flow, one may not be interested in the microscopic 
properties but merely in certain macroscopic properties. In the control 
problem, the stochastic model distinguishes between open and closed 
looped systems but the deterministic model does not. Another reason for 
using a stochastic model may be that this model is closer to the physics 
of the actual situation. 

In any case the purpose of this section is to lay the ground work 
for stochastic control theory. Our attention will be focused on certain 
concepts of stochastic processes and random differential equations. 

To describe a stochastic process rigorously would require measure 
theory and a great deal more time. Our approach will therefore not be 
rigorous, but hopefully it will be complete enough to get across the 
basic ideas. For the rigorous approach, see either Doon or Gikhman and 
Skoralnodses 

A real random variable € is a set of numbers or events together 
with a probability measure defined on this set. It is characterized by 


its distribution function F(x) which is defined by 


Ge) = Pte < ox} 


“Nee. J. L., "Stochastic Processes,'' Wiley, Inc., New York (1963). 


MU ecikdbaneermy, I. I. and A. V. Skorokhod, "Introduction to the Theory of 
Random Processes,'' W. B. Saunders Company, Philadelphia, Pa. (1969). 
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where P {& <x) is the probability that § is less than or equal to x. 
The distribution function is nonnegative, nondecreasing, and continuous 
from the left; also F(- ~) = O and F(™) = 1. 

Analogously, if € is an n-truple of random variables, its distri- 


bution function is a function of n real variables. 
Ruste = < 6.00 
F(x); Xo» > x) iP {E as Bae 3 E < x} 


and F is called a joint distribution function of the variables Ke The 
function F(x), Koorees x) is uniquely defined in n-dimensional Euclidian 
space ED is non-decreasing, and is continuous from the left with respect 


to each variable. Furthermore, 


F(x), Kore Kin 7 Ms Xe pooees 


and 


Cope pape te eae) r= gl) Ge 


where po) denotes the distribution function of the i-truple (Eporee> E5)- 
A random function or a stochastic process is a random variable &(t) 
which is a function of time. As time varies, &(t) describes the evolu- 
tion of the process. If a random process is recorded as it evolves, the 
recorded function €(*) describes only one of the many possible ways in 
which the process might have developed. The recorded function &(°) is 
called a sample function of the random process. For each fixed value of 
t, the quantity §(t) is a random variable. 
Whereas a random variable is characterized by a distribution function, 
a stochastic process is characterized by a set of joint distribution 
functions. Assume that it is possible to assign a probability distribution 


to the multidimensional random variable ee = (E(t,)> E(ty)s--es g(t )) 


for any n and arbitrary times t.. The distribution function 
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F(x, , Koaeees X35 tyoeees t) = {E(t,) S Kjposse5 aCe) Se 


is called the finite-dimensional distribution of the stochastic process 
E(t). For F to be a distribution, it must satisfy the following com- 


patibility conditions: 


F(x), Korres Xo CO, Ee eas tyorees t) = F(x), Koorees XS tyorees t) 
for i < n and 
F(X sees Xi Cporees t) = ape ure a8 ida abe Ea) 
where dpo---od, is an arbitrary permutation of the indicies 1, 2,..., n. 


The mean value of a stochastic process is defined by 


co 


m(t) = E[&(t)] = | x d F(x, t) 


—CO 


where E is the mathematical expected value. The mean value is thus a 
function of time. Higher moments of € are defined similarly. 


The covariance of the stochastic process is given by 


Ta (Sent) 


Covn Ge) ey ACS) 5 = SEINE) em) i CaiCs)) eer (Sy Bl 


ff (Gam Ce) es Gyan (S)) ancl ees Canyons) 


Our definition of a stochastic process is very general, and 
most systems which come under this definition would be mathematically 
unmanageable. Some specialization of the theory which makes it possible 
to characterize the distribution of E(t,), E(ty)r-++5 g(t.) in a simple 


way are particularly attractive. For instance, if the distribution of 
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E(t, )>+--8(t,) is identical to the distribution of E(t, i 
E(t, + T)y +++ 58(t, + tT) for all T and all arbitrary choices of the 


times t tis then the stochastic process €(t) is said to be stationary. 


Backs 
Ihe oe first and second moments E[&] and BLES] of the distributions 
are equal, then the process is weakly stationary. 

Our discussion of control systems has been limited to systems in 
which knowledge of the system at time t together with the governing 
equations suffices to describe its future evolution. Knowledge of the 
past when the present is given is superfluous relative to the future 
evolution of the system. The stochastic system analogy of this situation 
is the Markov property for random processes; these are stochastic process- 
es in which the past and future of the processes are conditionally 
independent. In order to define a Markov process, the conditional 
probability and the transition probabilities have to be defined. The 
conditional probability P(A|B) is the probability that A will occur if B 
has occurred. Given a sequence of times th < to Sooo to Sie, wns 
probability that €(t) < x if the sample function €(*) has already taken 
the values E(t,)s E(ty),+++,8(t_) is denoted by P(&(t) < x[E(ty),.-+5 


E(t_))- A stochastic process is said to be a Markov process if 
P(E(t) < x[E(t,),---» S(t,)) = PEC) < x/E(e,)) 


The transition probability distribution F(x, Ely s) is defined by F(x, 
tly, s) = P(&(t) < x/&(s) = y). If a stochastic process is a Markov 


process, its finite distribution functions are given by 
F(x,, Koseres KS typos 


F(x,5 t,) F (x55 ty |x,> t,)---FQ, Ea Eaare ew 


45 


This results from an application of the Baye rule. A Markov process is 
thus defined by two functions, the absolute probability distribution 
F(x, t) and the transition probabilities F(x, tly, s). 


Consider a system with the following dynamic equation: 
2 (ES > WY) or Ew (e) (4.1) 


where € is a small parameter and w is a stochastic process. Since w is 
stochastic, the state of the system x will also be stochastic; thus, we 
are interested in solving stochastic differential equations. Further- 
more, our interest is not with a particular sample function x(*) which 
is a particular discription of the state of the system during one run 
through the process; our interest is with the statistical properties of 
the stochastic process x (t). 


Consider the linear stochastic differential equation 
dx = A x dt + dw (4.2) 


where w is a stochastic process. In order to make some progress in 
finding the statistical properties of x, assume that w is a Wiener 
process. 

A Wiener process is a Markov process which satisfies the following 
conditions: 


1. It is a second order process; that is, for all t 
AN 
ale Ce) < © 
Hence, the mean m(t) exists as well as the covariance function 


r(s, t) = cov [w(t), w(s)] 
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2. The process has independent increments; that is, for arbitrary 
Ene Se a ee Eh emsinernements 
ib Z n 
x(t) - x(t 3), x(t,_j) - x(t,_5),---sx(t,) - x(t,), x(t) 


are independent.* 

3. The distribution of x(t) - x(s) for arbitrary t and s depends 
only on t-s. In this case, the process is said to have stationary 
increments. 

4. The transition probabilities are Gaussian. In the one- 


dimensional case, the transition probability density is 
oe ++ Ne, wiht, O) = E exp - we /2At 
4 7 V2TAt 


5. w(0) = O with probability one, and E[w(t)] = 0 for all t>0O. 

Sample functions of a Wiener process have interesting properties. 
They can be continuous functions but are nowhere differentiable. Their 
paths are of infinite length. Yet it is for just such perturbations 
that (4.2) will be solved. 

If w in (4.2) had bounded variation, the solution could be written 


in terms of the transport matrix (x, t) of the linear system 


y= Any, (4.3) 
The solution of (4.2) would be 
t 
Ed (i) sO CEO) eect | O(t, tT) d w(t) (4.4) 
0 


where the value of x at t = 0 is the random variable c. The expectation 


of c is m and its covariance matrix is [.. 


*Independent random variables are defined on page 7 of Masi.” 
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The integral 


t 
O®(t, T) d w(t) 
0 


is a stochastic integral. Since the transport matrix ®(t, T) is 
deterministic and has continuous derivatives, one way of defining 


this integral is through integration by parts. 


t 
| O(t, tT) d w(t) = O(t, t) w(t) - S(t, 0) wO) 
0 


t 
a® 
- J el Ce) awiGe) adar 


It follows from (1.15) and other properties of the transport matrix 
that 
t 


E 
| Oe, wt) a wie) = we) = ©Ces ©) w(O) a | O(t, tT) A(t) w(t) dt 
0 0 

(4.5) 


The integral on the right exist for almost all sample functions since 
the sample functions of w(t) are almost all continuous. This way of 
defining the integral has the desirable feature that the integral can 
be interpreted as an integral of sample functions. It does not, how- 
ever, preserve the intuitive idea that the integral is a limit of sums 
of independent random variables nor can it be extended to the case 
where ® is stochastic. Doob gives a more formal definition of the 
integral together with detailed proofs of its stochastic properties. 


The expected value of this integral is computed as follows: 
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t 
B, | f CC) ad vo | = E[w(t)] - ®(t, 0) E[w(0)] 
0 
va 
0 


t 
=m(t) - O(t, 0) m(O) + { OES) ANG) eam! Ge) acic 
0 


Ol Can ye nA Ge) aa 


Hence 


t t 
E i OlG=aes) mec vo] = | OES eaccleanmn Gr) (4.6) 
0 


The properties of the solution of the stochastic differential 
equation (4.4) will now be investigated. Since x is a linear function 
of a normal process, it is also normal and can be characterized com- 
pletely by the mean value function and the covariance function. Since 


the expected value of the Wiener process w(t) is zero, 


t 


O(t, 0) Elc] + E | O(t, T) d w(t) 
0 


E[x(t) ] 


O(t, 0) My 


where My is the expected value of the initial condition c. Hence 


m, (t) = E[x(t)] = @(t, 0) My (4.7) 
Taking derivatives yields 


aoe = GE PE, 0) my = ACE) (te, 0) my = ACE) m, (4.8) 
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Thus the mean value satisfies the linear differential Equation (4.3). 
The covariance matrix is more difficult to compute. In order to 
simplify the calculations, assume My = 0; hence, E[x(t)] = 0. This can 


always be achieved by subtracting mi, OMNES LOGS mame, 


Ris, eaicov Ix). x(t) = EleG) x (1 


Sy 
=i lee, t) x(t) + I d(s, 0) d “| x! (t) 
iE 
fs) 
= 6(s, t) Elx(t) x (t)] + | ®(s, 0) Eld w() x (t)] 
a 
=Ol (shat) eR (tema) (4.9) 


The integral is zero since w(0) and x(t) are independent for s > t. 
Set PG) — Rene) = Eee) x (t)]. Then P(t) is the variance and is 


therefore the function of interest. 


t 
E Ce 0) c+ | Oey Od a) 


0 


t alg 
oF (se O) cect { OGEs Oo) ad “) 
0 


ace, O) Ble &] © '@, 0) 


P(t) 


iC 
+ o(t, 0) E|c { d wi(a) o(t, 9) 
0 
t 
+ { ®(t, tT) Eld w(t) c | @(t, 0) 


(s) 


i 1E 


i | I ite, 2) Bld @@) da @)l Ot, o) 
OnaO 
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The increments of the Wiener process are independent of C; hence 
E [ec d’ w(G)] = E [d w(t) c1] = 0 
Moreover, from the properties of the Wiener process 
il 
E [d w(T) dw (o)] = 0 
if dt and do have no parts in common; otherwise 
ae 
E [d w(T) dw (t)] = Ro dt 


where Rk, is the covariance matrix of the Wiener process w. The final 


expression for P is then 


t 
P(e) = OC, O) F OI, OY | o(t, T) R(t) oi (t, t) dt 
0 


(4.10) 


A differential equation for P can be obtained from this expression 


for P simply by differentiating 


dp fa T di 
@. Ee a(t, 0)| r o"(t, 0) + o(t, 0) roe, 0) 
fe 
+ O(t, t) R(t) Ot, @) 2 il ee R(t) Oo (, 1) aL 


0 


t 
+ | G(r, 1) R(t) SO (t, 1) at 
0 


Sil 


The transport matrix satisfies 


a® A 
me (ent) = A O(E, @) 
and 
OG, ©) — Al T 
ae © (Ce, ©) A 
Hence 
ma =A O(t, 0) To (t, 0) + O(t, 0) TO (t, 0) A 


ie 
+ R(t) + J A O(t, T) R(t) oi (t, 1) dt 


(c 
A { O(t, T) R(t) oC. a) AD at 


iE 
O(t, tT) R(T) 6 (t, 1) af 


Ge h [ee 0) Fete, O) { 
0 


t 
-| O(t, 0) To (t, 0) + | O(E) 0) RON (ey al |! At 
0 


+ R(t) 


Thus from (4.10) 
dP _ T 
—— S A IP op 1 AND oP R (ania) 


P(O) =f (4.12) 
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THE KALMAN-BUCY FILTER 

The solution of the optimal control problem for a linear stochastic 
system is given by the separation theorem. It consists of an optimal 
filter for estimating the state of the system from the observed data and 
a linear feedback of the estimated state of the system; see Figure 7. 
The linear feedback is the same as the feedback that would be obtained 
if there were no stochastic perturbation of the system. This section 
will develop the explicit computational schemes for solving the filter- 
ing problem. 

Suppose we have the stochastic process described in the previous 


section 


dx = A x dt +d w(t) Gr BD) 


x(0) =e (5.2) 


where w(t) is a Wiener process and c is a Gaussion zero mean n-vector. 
In an actual case in which the process is realized, it is important to 
know the state of the system. It is, however, not always possible to 
measure x directly; instead, a set of quantities z(t) dependent on x are 
measured. Assume that the dependence of z on x is linear and is given 


by 


dz = H x dt + dv Ges) 


where the perturbation v is a Wiener process independent of x. 

The filter problem can be formulated as follows. Assume that a 
realization of the output z has been observed over the interval 
0 <T< t. Determine the best estimate of the value of the state vector 
x at time t. It is assumed here that the admissible estimates of x are 


linear functionals F(z) of the observed output z. The criterion 


5)3) 


for determining the best estimate is that the mean square estimation 
error be a minimum. This best estimate X(t) is dependent on the values 
of 2(@) yin gehe antervalyOR< sa < see ands elicansibeRprowed senatptyicma 


linear combination of the values of z on this interval. 


t 
R(t) = i Kes t) ol BC) (5.4) 
0 
Since z(t) is a stochastic variable, &(t) is a stochastic integral. 
Interpolation and extrapolation are two problems that are related 
to the filtering problem. The interpolation problem is one of estimating 
the state at some time T < t; the extrapolation problem is one of esti- 
mating it at some time tT > t. This latter problem is the one which is 
of interest to the stock market investor. 
The condition that X(t) is the best estimate from among all linear 
functionals of z(t) for the state vector x in the least squares sense is 
stated mathematically as follows. For every constant vector i and 


linear functional F, 
a D 
BL{ATGe(t) - (t))}7) < ELA (x(t) - F(2))}] — .5) 
where all variables have a zero mean. 


aISRCE) I) = BibaCe)] = BLE@) I = © 


Now set 


~ nw 
X=x- X 
where x is called the minimum error vector. 


we D7 <a Ee & @@) = SIA 


2 HO 37] & IA BA" GG) = DJ 


BO! @@) = D1 
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For all A and F(z), the criterion (5.5) requires 


me! Ge) = Syl = at = GG 2S S60 


This can be true only if 


o2nol 2. Ge =D eM aee@ = on 


But this implies that 


me @@) =a "jy <0 


for any linear combination F(z) of elements of z; hence 
ll 
E[x F (z)] = 0 (5.6) 
An integral equation for the kernel K(t, T) can be derived from 
(5.6). This kernel is not a stochastic quantity, and it can be de- 


termined independent of the realization z(*). For F(z) = z(t) - 


ZO) Ee ONO sete ee expressions (5)16)) sy lelds 


Tin) @@) =2@1 = nee) G@) = 2) 71 


e T 
E | X(t) || H(s) x(s) ds +d vo} 


O 


fe TG T 
i Ke, 7) 120 | {J H(s) x(s) ds + wo 


ll 
ea) 
———_ 


0) Oo 


t T 


te 40 
=6 | | Cy GGG ae deci (H(G)) <(G)de cera())in 
0 


2 | ee, 1) HG) ee) = Gy EE@) ds oe 
0 % 
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ie iG 


a i | Kes HG) ss) dv (s) dr 
0 


uF { Kes) idviGs) x! (s) H. (s) ds 
0 Oo 


t TE 
+ | | K(t, r) dv(r) ee) 
0 fe) 


From the properties of Wiener processes, 


i iG E 


E { | K(t, r) dv(r) aa (a) = i K(t, s) R (s) ds 
oO we Oo Mi 


where R is the covariance matrix of the process. Furthermore, 


dv(s) and x(s) are independent, so E dv(r) ae) = 0; hence 


ie 


as 
Elx(t) (2(t) - 2(0))"] = | | | K(t, x) H(r) Elx(e) x (s)] 
fo) 0 
(5.7) 


Et (e) dr + K(t, s) ol ds 


for all o and t. On the other hand, from (5.3) 


E T 
Bint) @@) = B@))"] =m Ee { ils) «(e) ds + vio | 


0 


T 


| BinGe) <a) ECs) ae 
oO 


i is 
= i i K(t, r) H(r) Elx(r) x'(s)] H'(s) dr 
0 


oO 


+ K(t, s) RCs) ds 
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where the last equality results from (5.7). Since this equation holds 


Eom allluronande ct sinyteheysimitcersvalen[|Onnatells 


i 
1 T , a T 
Kes Ss) R(s) = E[x(t) x (s)] H (s) - | K(t, r) H(t) E[x(r) x (s)] H (s) dr 
0 
(5.8) 
This is a nonhomogeneous integral equation for K(t, s). Its kernel 


is H(r) E[x(r) i) E-(@)- Since it corresponds to a positive definite 

quadratic form, all its eigenvalues are positive and the equation has a 

solution. Unfortunately, it is not possible to calculate K(t, s) from 

this equation because E[x(r) <(@)1)- the covariance of x(s), is unknown. 
A different equation for K(t, s) can be obtained from (5.8) by 


differentiating both sides of it with respect to t. 


oK(ES) p(s) = S- Elx(t) x"(s)] Hs) 
= K(e, t) H(t) Elx(t) x! (s)] H (s) 
t 
{ o Reet) HG) BERG) ©] &G) de 
0 t 
By (5.1) 


dx = Ax dt + dw 


Hence 


Diente) (ey) HUG) = AG) ie) 2G] ae ae 


+ E[dw(t) <i 
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The second term vanishes since dw(t) and es) are independent if 


ESL hel Sieyeicenicis 


cae Rs) = [A(e) — Ke ©) (ey) Elx(e) x'(s)] H'(s) 


ie 


es | a Ge r) H(r) E[x(r) (Ss) I H! (s) dr 
0 t 


Use of the integral equation (5.8) to obtain an expression for 


E[x(t) x e)j H (s) yields 
[2K 2) + eee, 0) HCE) KOE, 8) — ACE) KEE, 99] RYO 


ie 
ai [28g + K(t, t) H(t) K(t, r) - A(t) K(t, o| 


0 
H(x) Elx(r) x (s)] H’(s) dr (5.9) 
Set 
Die, s) = {EE 8) 4 Kee, ty H(t) K(t, s) - A(t) K(t, 5) 
Then by (5.9) 
t 
Y(t, s) = - I Ue, x) Wr) Elx(e) x"(s)] W'(s) RUv(s) dr 
0 
(5.10) 


Since the kernel of this integral equation corresponds to a positive 


definite quadratic form, the only solution of (5.10) is 


UGS) =O 
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This yields the following differential equation for K(t, s): 


: ae S) = a(t) K(t, s) - K(t, t) W(t) K(t, s) (5-12) 


From the integral equation (5.8) for K(t, s) 


K(t, t) R(t) = Elx(t) x (t)] H(t) 


ie 
- I K(t, 1) H(t) Elx(x) x! (t)] wi (t) dr 
0 


On the other hand, 


t 
B[&(t) x! (t)] = E ! K(t, x) dz(r) e] 


0) 


t 
= | KG) HG) Bilas Ge) x! (t)] dr 
0) 


t 


+ | K(t, r) E[dv(r) eI 
0 


where the second integral vanishes; hence 


K(e, £) R(t) = Elx(e) x (e)] W(t) - ELR(t) x (t)] WC) 


El(x(t) - (t)) x! (t)] H(t) 


= {E[X(t) &1(t)] + ELR(t) #2 (t)]} w(t) 
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The condition that ¥ is the best estimate from among all linear func- 


tionals of z(t) leads to the result that E[x(t) aH ayy = 0. Hence 


P(e) EXC) S Ble) =2¢z)] wom) = R(t, &) Ra GE) 


From the stochastic integral, 


t 
(CE) = i K(t, r) dz(r) 


0 
iE 
d&(t) = K(t, t) dz(t) + | SRE) ae (r) at 
0 
(e 
=P HR) a(t) + | (A(t) K(t, x) - K(t, t) H(t) K(t, r)) dz(r) dt 
0 
or 


Ga(e) = Ate) SGe) de eo Ae Rot (d2(t) = Ge) Ae) ae) G13) 


Since z(t) and presumably dz(t) are known, this is a stochastic differ- 
ential equation for X(t). 


Note that 


dz(t) - H(t) &(t) dt 


dz(t) - H(t) x(t) dt + H(t) x(t) dt 


dv(t) + H(t) x(t) dt 


From this expression and (5.13), we get the following stochastic 


differential equation for x 


60 


dx(t) - dx(t) 


a 
bald 
i 


=i] ES 
Axdt+dw-Adt-PHR (dv + H(t) x) 


at till 


Ni de 4 dw — PH Ras dv - PH RL H x dt 


[A= pe a a] & des@n ep ae ee ae 


with x(0) = x(0). By the methods developed in the previous section 


for stochastic differential equations, 


P(t) = E[x(t) x2 (t)] 


ie 
a(t, 0) To (t, 0) + { (t, 0) [Q(o) 
0 


+ 


P(o) HE (0) R-* (0) HG) Bo] OE. o) de Gil) 


where ®(t, T) is the transport matrix associated with the linear 


differential equation 


dy a a i il! 
ae (A - PH RL H) y Ges) 
and where 
Q dc = Bide) Gy Gy (5.16) 
from (4.11) 
dP i 1 T 


SS = Pek ae H) P + P(A - P He Ro H) 


£ Oe) put Se H Pl 
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or 


dP IL 


oT AP+PA-P HR 


ae Bl BP = @ Gal 


EGO) eal (5.18) 


This set of equations finishes the solution of the filter problem. 

The optimal filter is a feedback system which is described by the 
stochastic differential equation (5.13). It is obtained by taking the 
measurements z(t), forming the error signal z(t) - H(t) &(t), and feed-— 
ing the error forward with a gain P(t) H(t) mC) P(t), the error 
variance, is obtained as a solution to the nonlinear Riccati-type 
equation (5.17), H(t) is a known transformation matrix, and RY is 
the variance of the Wiener process dv. A block diagram of the filter is 
shown in Figure 8. The variables appearing in this diagram are vectors, 
and the boxes represent matrices operating on vectors. The double lines 
which indicate direction of signal flow serve as a reminder that multiple 


signals rather than a single one are being directed. 


Figure 8 -- Optimal Filter 
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